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Using the Reduced Basis approach, we efficiently compress and accurately represent the space 
of waveforms for non-precessing binary black hole inspirals, which constitutes a four dimensional 
parameter space (two masses, two spin magnitudes). Compared to the non-spinning case, we find 
that only a marginal increase in the (already relatively small) number of reduced basis elements is 
required to represent any non-precessing waveform to nearly numerical round-off precision. Most 
parameters selected by the algorithm are near the boundary of the parameter space, leaving the 
bulk of its volume sparse. Our results suggest that the full eight dimensional space (two masses, two 
spin magnitudes, four spin orientation angles on the unit sphere) may be highly compressible and 
represented with very high accuracy by a remarkably small number of waveforms, thus providing 
some hope that the number of numerical relativity simulations of binary black hole coalescences 
needed to represent the entire space of configurations is not intractable. Finally, we find that the 
distribution of selected parameters is robust to different choices of seed values starting the algorithm, 
a property which should be useful for indicating parameters for numerical relativity simulations of 
binary black holes. In particular, we find that the mass ratios mi/mi of non-spinning binaries 
selected by the algorithm are mostly in the interval [1,3] and that the median of the distribution 
follows a power-law behavior ~ (rai/m2)~ 5 ' 25 . 



I. INTRODUCTION 

The upcoming generation of advanced-sensitivity 
ground-based gravitational wave interferometer detectors 
(i.e., advanced LIGO, advanced Virgo, Indigo, and KA- 
GRA) [IH1] brings an increasing demand to accurately 
and efficiently represent gravitational waveforms from 
generic precessing compact binary sources [5HT2"]. Such 
waveforms for quasi-circular inspirals are parameterized 
by a set of eight intrinsic physical quantities - two masses, 
two spin magnitudes, and four spin orientation angles on 
the unit sphere (8D) [32] . 

The large dimensionality of this parameter space 
makes gravitational wave searches, parameter estimation, 
and modeling prohibitively expensive and computation- 
ally unfeasible with most methods. This problem is called 
the "curse of dimensionality" [13] and, in particular, is 
a major hurdle for modeling astrophysical gravitational 
wave sources. In this setting, the cost of numerical rela- 
tivity simulations of the full Einstein equations describing 
the inspiral, merger, and ringdown of binary black hole 
coalescences is so expensive that an optimal or nearly 
optimal criterion for selecting which points in parameter 
space to simulate is thus critical. This is also the case 
when numerical simulations are used to calibrate or build 
effective one body [T1HT7] or phenomenological [TM2"U] 
models. 

To address these issues, the Reduced Basis (RB) ap- 
proach was introduced to gravitational wave physics in 
Refs. [2TJ [22] and shown to efficiently compress the space 
of waveforms with a very small loss of accuracy, typically 
of the order of numerical round-off, for any given bounded 



parameter domain. The compression is accomplished by 
determining a set of nearly optimal physical parameter 
values from which a basis is constructed to represent any 
given waveform within the same physical model through 
its projection onto this basis. Details of the algorithm, 
which is known as the greedy algorithm can be found in 

[Sung. 

The Reduced Basis approach has several uniquely ap- 
pealing features such as its hierarchical compression of 
the waveform space, its ability to handle large numbers 
of parameters, and its ability to identify the most relevant 
points in the parameter space. The latter may be par- 
ticularly important for guiding which numerical relativity 
simulations to perform in order to calibrate effective-one- 
body or to fit phenomenological models. 

In this paper we apply the RB approach to non- 
precessing binary black hole inspirals that are non- 
spinning (2D), have equal spin magnitudes (3D), or have 
unequal spin magnitudes (4D). To our knowledge, this 
is the first paper, along with Refs. [35] [23], applying re- 
duced order modeling techniques to gravitational waves 
with a larger number of parameters than previously con- 
sidered. 

The results presented in this paper point towards the 
possibility that the eight-parameter space of waveforms 
for precessing binary inspirals might admit a remarkably 
compact representation. For the non-precessing inspi- 
ral waveforms, we find a very moderate increase in the 
number of RB elements as the dimensionality of the pa- 
rameter space is increased from 2D to 4D. While this 
number could and might change significantly for precess- 
ing binaries, former results in |24) show dimensionality 
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reduction in the precessing case (with respect to the spin 
dynamics). 

II. SPIN-ALIGNED POST-NEWTONIAN 
WAVEFORMS 

Throughout this paper we use the so-called restricted 
TaylorF2 post-Newtonian (PN) waveforms (see [25] and 
references therein for a more detailed discussion of Tay- 
lorF2 and other PN approximants). These waveforms use 
the stationary-phase approximation to construct analytic 
frequency-domain waveforms of the form 

h(f) = Ar 7/ *e i *W ) (1) 

where *f?(f) is expressed as a polynomial in the PN ex- 
pansion parameter v — (ttM/) 1 / 3 and A is a constant 
such that the h is normalized to unity. Spin-independent 
corrections to ^(/) are currently known to v 7 , or 3.5PN 
order [26 . The dimensionless spin vectors of the com- 
ponent masses, Xi = Si/mf with < \xi\ < 1, can also 
enter the corrections to ^(f) and these spin-dependent 
corrections are known to 2.5PN order. These include 
spin-orbit corrections of the form Xi ' Lpf (where Ljy is 
the "Newtonian" orbital angular momentum, which is 
transverse to the orbital plane) at 1.5PN and 2.5PN or- 
der, and spin-spin and self-spin corrections of the forms 
Xi • X2 and xf at 2PN order [27H5U]. 

For a generic binary system the spins and orbital plane 
precess, so the relative orientations among x\, xi an d 
L N , and thus the spin-dependent PN corrections of &(f), 
vary on a precessional timescale. Currently, there are no 
frequency-domain template families for precessing spins, 
although there is at least one such effort underway [5T] . 
In the "spin-aligned" case, when each spin is aligned or 
anti-aligned with Ln, the binary does not precess and 
the relative orientations of xi, Xn an d Ln (and the spin- 
dependent corrections to &(f)) remain constant through- 
out inspiral. Thus the simple TaylorF2 model is valid and 
does not need to be augmented with precession equations. 
The spin-dependent corrections are expressed in terms of 
the (constant) projection of each spin along the orbital 
angular momentum: — 1 < Xi = Xi ' < 1- 

In the restricted TaylorF2 approximation, the ampli- 
tude is expanded only to leading order while the phase 
is expanded to a higher order. In all cases, the spin- 
independent contributions to the phase are here included 
up to 3.5PN order, as given by Eq. (3.18) of [25J- In cases 
where spin is included, the 3.5PN non-spinning phase is 
augmented with spin-orbit, spin-spin and self-spin cor- 
rections through 2.5PN order, as can be read from Eqs. 
(6.22)-(6.25) of [3D]. 

III. RESULTS 

We consider binary black holes [33] with individual 
component masses rm € [3, 30] M© and dimensionless 



spin magnitudes Xi € [ — 1,1]. Each mass direction (or 
dimension) in the training space (any discretization of 
the continuum space of parameters) is sampled with n m 
points and each dimensionless spin magnitude direction 
with n s points. After some numerical experimentation, 
we found that the majority of the selected mass com- 
ponents have values at the lower end of the considered 
range. Thus, we use a training space with points loga- 
rithmically spaced in the mv-m% plane to provide a suf- 
ficiently dense cluster of points in the training space at 
those lower mass values. The inner products between 
any two waveforms are weighted by the reciprocal of the 
power spectral density (PSD) for advanced LIGO [33] as 
given by the fitting formula in Ref . [32j , with a lower cut- 
off frequency of 10Hz and a maximum one correspond- 
ing to the highest ISCO frequency. For clarity, when 
we refer to the number of RB elements we mean the 
number of basis vectors needed to represent the training 
space to within the specified tolerance, typically round- 
off (~ 1CT 14 - 1(T 12 ). The number of RB elements also 
equals the number of selected parameter tuples of mass 
and spin magnitudes. 

Compared to the non-spinning case, the extra dimen- 
sionality of the spin parameter space (for non-precessing 
inspiral waveforms, as considered here) has a remarkably 
small effect on the final number of RB elements needed to 
represent the entire training space with numerical round- 
off precision. In fact, as will be discussed in Section |lII C[ 
the same RB represents any waveform, and not necessar- 
ily a member of the training space, in the same range 
of parameters to within essentially round-off precision. 
These are the most important results of this paper. 



A. Effect of increasing dimensionality 

Figure [I] shows the maximum representation error for 
training spaces, corresponding to 2D, 3D, and 4D pa- 
rameter spaces of non-precessing inspiral waveforms, as 
a function of the number of RB elements. For each wave- 
form h in the training space we compute the representa- 
tion error [45] \\h — P/v[/i]|| 2 and then plot its maximum 
over the entire training space, as the number of RB el- 
ements is increased. In all cases Fig. [T] shows a similar 
behavior, namely, an initially slow fall-off in the represen- 
tation error and a later rapid convergence to round-off. 

The 2D runs were done with up to n m = 400 points in 
each mass dimension, for a maximum size of the training 
space of = 1.6 x 10 5 total samples. For the 3D (4D) 
runs we used up to n m — 200 (100) and n s — 50 (20) 
samples in the spin direction for a training space with 
a maximum size of 2 x 10 6 (4 x 10 6 ) samples. In order 
to directly compare the 3D and 4D cases in Fig. [TJ we 
use n m = 100 samples for both since, as discussed in Sec- 
tion [lire] more samples are not needed for high accuracy 
reconstruction of arbitrary 3D or 4D waveforms. 

The key result from Fig. [I] is the very moderate in- 
crease in the maximum number of RB elements needed 
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as the dimensionality increases from 2D to 3D to 4D. At 
an error of 10 -11 the 2D case requires « 1,725, the 3D 
case « 1,824, and the 4D case s» 1,839 basis vectors. The 
number of basis vectors that result by increasing the di- 
mensionality of the parameter space from 2D to 4D is 
only a 6.6% increase and from 3D to 4D requires only 
15 more basis vectors. Despite doubling the dimension- 
ality of the parameter space, the number of basis vectors 
needed to represent the four-parameter waveform space 
is nearly the same as for non-spinning binaries. 

Naive scaling arguments would suggest an increase in 
the number of RB elements by orders of magnitude. If the 
2D case requires ~2x 10 3 basis vectors then one might 
be lead to think that doubling the number of parameters 
could require more by a factor the size of the training 
space in the spin directions, which would here be about 
20 2 for a total of ~ 8 x 10 5 . The fact that we observe 
an exceedingly small relative increase in the number of 
RB elements lends hope that the number of basis vectors 
needed to represent the full 8D parameter space for inspi- 
ral waveforms is also much smaller than what might be 
anticipated by estimating the volume of the parameter 
space or by using naive (such as equally spaced or ran- 
dom) sampling techniques. This is only suggestive since 
precession has significant effects on the structure of the 
inspiral waveforms 6, 27, 30 . But if it holds, the curse of 
dimensionality might be beaten. In fact, indications that 
the problem, in the presence of precession, is amenable to 
dimensional reduction, have already been found through 
a Principal Component Analysis of the precessing dy- 
namics of compact binary inspirals 24 . There it was 
found that for the case of a random selection of precessing 
binaries with the same total mass there are three combi- 
nations of spin orientations that are semi-conserved (in a 
statistical sense) throughout the inspiral. The presence 
of such combinations implies that the dynamical config- 
urations can be parameterized by a reduced number of 
independent parameters than the fiducial seven (1 mass, 
6 spin components). 



B. Parameter values selected 

The Reduced Basis-Greedy Algorithm approach pos- 
sesses several unique features. One of them is being able 
to identify, in a precise mathematical sense , a nearly op- 
timal set of points in the physical parameter space |33j . 
In this section we analyze the structure of the selected 
parameter values for the models and scenarios here stud- 
ied. 



1. Mass parameters (2D, 3D, 4D) 

Figure [2] shows the distribution of selected points for 
the mass ratio mxjmi in the non-precessing (4D) case. 
There we plot log 10 (7711/7712) to obtain a symmetric rep- 
resentation in mi and 7772 . Since in this paper we focus on 
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FIG. 1: The reduced basis representation error as a function 
of the number of reduced basis for different dimensionalities. 
The greedy error shows the rapid exponential convergence to 
round-off. The number of reduced basis elements increases 
very mildly from the 2D (no-spin) to the 3D (equal-spin) case 
and to the 4D (generic non-precessing spins) case. This shows 
that the increasing dimensionality may not be a major obsta- 
cle for the Reduced Basis approach. 



the range mi, mi € [3, 30] M , the mass-ratio has a range 
777,1/777,2 € [1/10, 10]. The selected values by the greedy al- 
gorithm have a strong peak at the equal-mass case. This 
follows from the selected individual mass components, 
which are clustered at low values (see also Fig. [3]) . 
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FIG. 2: Histogram of the selected parameter choices for the 
mass-ratio m\Jm^ in the 4D case. Note the large spike near 
equal masses and the rapid fall-off to small counts for increas- 
ing disparity in the component masses. 

We have found similar results for the distribution of 
mi/777,2 in the 2D and 3D cases, also resulting from the 
low-mass selection. This is also visible from Fig. [3j which 
shows the selected parameter values in the 777,1-777,2 plane. 
Notice that there is little qualitative difference in going 
from 2D to 3D to 4D in Fig. [3] This indicates that the 
problem might admit some kind of dimensional separa- 
bility. For example, it might be possible to initially select 
(7711,7712) values through a 2D greedy strategy and then 
use only those to populate the m\-m2 plane in the con- 
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struction of very efficient and compact training spaces 
for building reduced basis in higher dimensional param- 
eter cases. This is beyond the scope of this paper, but 
it might provide a dramatic reduction in training space 
sizes and, if so, a possible avenue for expansions to higher 
dimensions and the full 8D problem. 
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FIG. 3: The reduced basis parameter choices in the m\-mi 
plane. The overall structure is similar for all three (2D, 3D, 
4D) non-precessing cases. The grey-scale bar utilizes a loga- 
rithmic scaling. 



2. Mass and spin parameters (3D) 

The 3D equal-spin parameter space (mi, ni2, Xi — X — 
X2) can be visualized directly. In Fig. Hwe plot the pa- 
rameter values selected for the reduced basis. Comparing 
this to Fig.[3]one sees again the usual selection of primar- 
ily low-mass systems. For the lowest mass configurations, 
a large range of spin values is selected. The figure also 
shows that very little of the bulk of the parameter space 
volume is used. The greedy algorithm mostly picks val- 
ues at some edges and corners in the mx-ra^ plane and 
preferentially picks binaries with both spins anti-aligned 
to the orbital angular momentum (x = —1); see also 
Fig. 5(a) showing a histogram of the selected equal-spin 



magnitudes. 

The global aspect of the RB greedy algorithm allows 
an identification of the underlying sparsity clearly visi- 
ble in Fig. [2J This provides another hint that: 1) the 
binary coalescence problem is amenable to dramatic di- 
mensional reduction with respect to the physical param- 
eters describing it, and 2) specific techniques for train- 
ing space construction exploiting sparsity might provide 
dramatic computational cost advantages when applied to 
higher dimensional problems (we elaborate more on this 
in Sec.lrvl). 



3. Spin parameters (3D, 4D) 



Fig. |5(b)| shows the 1,839 spin magnitudes, xi an d 
Xii selected by the greedy algorithm in the generic non- 
precessing 4D case. The training space corresponding 
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FIG. 4: The reduced basis parameter choices for the 3D case 
(rai,m2,x)- Comparing to Fig. [3] one can also see the se- 
lection of primarily low-mass systems. For the lowest mass 
systems a large number of spin values are selected. Few sys- 
tems from the bulk volume of the parameter space are chosen. 



to this figure has n s — 20 points in each spin direc- 
tion. As found in the 3D case, binaries with both spins 
mostly anti-aligned with the orbital angular momentum 
are more relevant than those that are aligned. 




(a) Spin magnitude 
distribution (3D case) 



(b) Spin magnitudes 
distribution (4D case) 



FIG. 5: Fig. (a) shows a histogram of spin values in the 3D 
case. Spins anti-aligned with the orbital angular momentum 
are predominantly chosen. Fig. (b) shows a two-dimensional 



histogram of spin values in the 4D case. Similarly to Fig. |(a)| 
anti-aligned systems \i = — 1 are predominantly chosen. The 
grey-scale bar utilizes a logarithmic scaling. 



The figure is asymmetric in xi 1 Xi f° r t ne following rea- 
son. Although the representation error of any two wave- 
forms is symmetric in (m\,Xi) ( m 2,X2), the greedy 
algorithm needs to select one set of parameters from the 
two possible choices introducing an asymmetry. This ef- 
fect is visually enhanced by the logarithmic scale used in 
the figure. 



As discussed in Sec. |IIIB 1| the distribution of masses 
selected does not change strongly as one goes from 2D 
to 4D. Similarly here, the distribution of spin values se- 
lected seems qualitatively similar in the 3D and 4D cases, 
suggesting that it might be possible to construct training 
spaces in a more efficient manner. 
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C. Monte-Carlo representation error studies 

To test the accuracy of the RB to represent waveforms 
that are not necessarily members of the training space, 
we randomly sample the corresponding parameter space 
used to build each basis, and for each sampled wave- 
form we compute its representation error by projecting 
the waveform onto the RB. We performed this test for 
all models (2D, 3D, 4D), but for defmiteness we show 
the more interesting 4D case in Fig. [6] Shown there 
is the computed waveform representation error for 10 7 
randomly chosen values in the 4D parameter space. The 
histogram shows that the overwhelming majority of sam- 
pled waveforms have a representation error near double- 
precision round-off values, with a few isolated cases with 
an error slightly larger than 10 -11 . 

For the particular reconstruction test shown in Fig. [6] 
we used a RB built out of a training space with n m 
1 00 points for each mass component and n s = 10 for 
each spin magnitude, as described in Sec. |III A| For the 
2D case, this number of mass samples is not enough to 
achieve such a low representation error. This indicates 
that there is a trade-off between sampling density in the 
mass and spin dimensions that can be used to reduce the 
representation error. The effectiveness of the reduced 
basis approach is such that, for example, even though 
a relatively small number of spin values are used in the 
training space to build the RB, a waveform with any spin 
value is represented within essentially machine precision. 




FIG. 6: Monte-Carlo Reconstruction error study for the 4D 
case. The histogram shows the distribution of the represen- 
tation error a = || h — P^[h]\\ for 10 7 randomly selected pa- 
rameter values. Notice the logarithmic scale on the vertical 
axis. 



D. Seed dependence 

Lastly, we investigate the dependence of the param- 
eters selected by the greedy algorithm on the choice of 



initial seed for the 2D case (non-spinning). The seeds 
are chosen to correspond to the 200 equal-mass binaries 
in the training space. For a given seed and tolerance er- 
ror, the greedy algorithm picks a set of nearly optimal 
parameters, or a greedy chain. Different seeds give rise 
to different greedy chains and the chosen parameters will 
not necessarily be part of other greedy chains. 

We showed in Ref. [55] that the representation error 
is robust to different choices of seed values by running 
the greedy algorithm multiple times, once for each pos- 
sible value in the training space for the seed [46!. Here, 
we discuss how the distribution of points is affected by 
different seed choices. 
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FIG. 7: Parameters chosen by the greedy algorithm (mi, left 
and mi/m.2, right) when using 200 seed values correspond- 
ing to equal mass, non-spinning binaries with different total 
masses. The lines are quantile curves representing the frac- 
tion of selected parameters with values greater than the value 
of the quantile. Black corresponds to the median and the 
bounding black curves correspond to the 0% (top) and 100% 
(bottom) quantiles. The red line is a fit to the median. De- 
spite the different choice of greedy parameters that are se- 
lected, the distribution of points is rather robust to the choice 
of seed. 

Fig. [7J shows the probability density over all 200 seeds 
for a given value of mi (left) and for mi/m.2 (right), 
where the mass ratio is taken to be such that if mi > m 2 
then the 1 and 2 labels are interchanged so that mi/m2 > 
1. The gray- level shaded lines correspond to quantile 
curves, with black indicating the median of the distribu- 
tion. Both plots indicate that the distribution of points 
is rather robust to different seed choices. This reflects 
the global nature of the greedy algorithm to select pa- 
rameters that are distributed similarly even though the 
selected parameter values themselves exhibit variations 
from one chain to another. These distributions also show 
that the majority of points selected, and overwhelmingly 
so, are those with low masses independently of the seed 
value. A similar distribution of points exists for the se- 
lected values of m?. The red line in the right panel of 
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Fig. [7] is a power-law fit to the median and is given by 
~ 1.83 (mi/m2)~ 5 ' 25 + 0.05. Most mass ratios selected 
lie within the range TOi/m 2 € [1,3]. 



IV. CONCLUDING REMARKS 

The main focus of this work has been the study of 
reduced basis in higher dimensional parameter spaces; 
more specifically, the study of non-precessing, spinning 
binaries. We have shown (see Figs, [j] and |6| that the 
number of reduced basis vectors needed to represent the 
full space of non-precessing inspiral waveforms with es- 
sentially round-off precision increases very mildly with 
the number of parameters from 2D to 4D. 

While the work presented here is limited to the non- 
precessing case, our results provide the first hints, to- 
gether with |24| for precessing inspiral dynamics (see 
also [H]), [H] for quasi-normal mode ringing, and stud- 
ies using a Singular Value Decomposition (SVD) tech- 
nique [23 EH EH] that a full representation of the eight- 
parameter space of inspiral-merger-ringdown waveforms 
(either numerical, Effective-One-Body, or phenomenolog- 
ical ones) might actually be achievable with a relatively 
compact reduced basis. If this is the case, then there 
is hope that the number of simulations needed to rep- 
resent the space of precessing binary inspirals might be 
relatively small (perhaps, on the order of several thou- 
sand for advanced LIGO but not orders of magnitude 
larger than this), thus allowing for a tractable number of 
numerical simulations if the parameters are chosen with 
malice of forethought from our reduced basis studies. 

Recently, SVD of gravitational waves [34] has been im- 
plemented in a realistic, low-latency search pipeline |36j . 



Such a pipeline could be used to improve upon the work 
in [37] to allow prompt follow-up searches for electro- 
magnetic counterparts to candidate gravitational wave 
signals in the advanced detector era. As a reconstruction 
and compression technique, RB has a number of similar- 
ities with SVD (e.g., output includes basis vectors and 
reconstruction coefficients in both schemes). Thus the 
RB results of this paper could be implemented in such 
a pipeline, in principle. Such practical data analysis im- 
plementations of RB will be left for future work. 

The exploration of the 8D parameter space will require 
further development and/or implementation of technical 
but critical aspects, including the efficient and adaptive 
sampling techniques for large training spaces (see, e.g., 
[38) ) and rapid evaluation of high accuracy quadratures 
for parametrized problems |39j . For the problem at hand, 



a splitting of dimensions as discussed in Sec. MB 1 might 
be useful. 
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